Proximity effects at ferromagnet-superconductor interfaces 
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We study proximity effects at ferromagnet superconductor interfaces by self-consistent numerical 
solution of the Bogoliubov-de Gennes equations for the continuum, without any approximations. 
Our procedures allow us to study systems with long superconducting coherence lengths. We obtain 
results for the pair potential, the pair amplitude, and the local density of states. We use these 
results to extract the relevant proximity lengths. We find that the superconducting correlations in 
the ferromagnet exhibit a damped oscillatory behavior that is refiected in both the pair amplitude 
and the local density of states. The characteristic length scale of these oscillations is approximately 
inversely proportional to the exchange field, and is independent of the superconducting coherence 
length in the range studied. We find the superconducting coherence length to be nearly independent 
of the ferromagnetic polarization. 
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I. INTRODUCTION 

In recent years, technological advances in materials 
growth and fabrication techniques have made it possible 
to create heterostructures including high quality ferro- 
magnet/superconductor (F/S) interfaces. These systems 
have great intrinsic scientific importance and potential 
device applications, includingpiquantum computers and 
magnetic information storageErQ. This has led to re- 
newed interest in proximity effects involving magnetic 
and superconducting compounds. Understanding how 
proximity effects modify electronic properties near F/S 
interfaces is constantly becoming more important as the 
rapid growth of nanofabrication technology continues. 

The juxtapositian of a ferromagnet and a supercon- 
ductor can results in a spatial variation of magnetic 
and superconducting correlations in both materials. The 
leakage of superconducting correlations into the non- 
superconducting material is an example of the supercon- 
ducting proximity effect. Similarly, the spin polarization 
may extend into the superconductor and modify its prop- 
erties, creating a magnetic proximity effect. 

In general, if one is interested in a microscopic solution 
of the F/S proximity effect problem valid at all length 
scales, oiie must solve the appropriale equations, e.g. the 
Gor'kovB, or Bogoliubov-de Genncsu (BdG) equations in 
a self-consistent manner and with as few approximations 
as possible. In practice, approximations are often made 
in the basic equations. Further, in many cases a sim- 
ple form for the pair potential A(r) is assumed, usu- 
ally a constant in the superconductor region, and zero 
elsewhere is used. Such crude non self-consistent treat- 
ments have been widely applied because of their simplic- 
ity. However they are valid typically only for length scales 
much longer than the superconducting coherence length, 
which characterizes the depletion of the pair potential 
in the superconductor near the interface, or in the case 



where the non-superconductor is very thinJH The super- 
conducting proximity-jcffect is linked to the phenomenon 
of Andreev reflectionEl. This is the process where at the 
interface, an electron is reflected as a hole, transmitting a 
Cooper pair into the superconductor and vice versa. The 
inhomogeneity in A(r) creates a potential well for quasi- 
particles, causing electon-hole scattering, and subsequent 
bound states below the maximum value of A(r). 

There are several quantities that can be studied, the- 
oretically or experimentally, in the context of character- 
izing proximity effects. The traditional descriptionQ of 
the superconducting proximity effect is through a char- 
acteristic proximity length which can be associated with 
the behavior of the pair amplitude^ F{r), the probabil- 
ity amplitude to find a Cooper pair at point r. This 
quantity does not vanish identically inside the non- 
superconductor. This is in contrast to the pair poten- 
tial A(r), which is of limited use, since it is zero inside 
the non-superconducting material unless it is arbitrarily 
assumed that a small nonvanishing pairing interaction 
exists there. An additional important quantity, which is 
now experimentally accessible thanks to improved STM 
technolog}^ which allow local spectroscopy to be per- 
formed, is the local density of states (DOS). This quan- 
tity reflects the one-particle energy spectrum, and there- 
fore one aspect of the proximity effect. 

For a nonmagnetic normal metal in contact with a su- 
perconductor, the proximity effect has been much studied 
and well understood for many yearsfl For clean systems, 
if the non self-consistent step function for the pair poten- 
tial is used, solutions lAitbe microscopic equations are rel- 
atively easy to obtain.E3"E^ Other approaches involve first 
simplifying the basic equations. One can, for example, in- 
tegrate out the energy variable in the Gor'koK equations. 
The resultant (quasiclassical) EilenbergeJli equations 
have the advantage of being first order, and therefore 
easier to solve. They can be extended to systems of ar- 



bitrary impurity concentration.lla Results that calculate 
the pair potential self-consistently are more sparse. TUa 
Eilenberger equations have been solved numerically,t2l 
and the DOS was calculated, with comparisons made be- 
tween self-consistent and non-self-consistent results. For 
systems in which the electron mean free path is much 
shorter than the superconducting coherence length, when 
the Eilenberger equaiions can be reduced to the—sim- 
pler Usadel equationalJ, a calculation of the DOSta has 
been performed. Numerical approaches which do not re- 
quire simplifying the starting equations are possible, al- 
though rare. Numerical self consistent solutions of the 
full Gor'kov ea^ations in heterogeneous systems have 
been obtainedJ13'E2l and from these the density of states 
and pair potential of normal metal-stmerconductor bilay- 
ers and multilayers were calculated. E3 

When the normal metal is replaced by a ferromag- 
net, the theoretical situation is much less satisfactory. 
The presence of the exchange field in the ferromagnet 
makes the overall physical and mathematical picture 
of the proximity effect in F/S systems quite different 
from its non-magnetic counterpart. Since on the mag- 
netic side Fermi surface quasiparticles with different spins 
have different wavevectors, numerical solution becomes 
much more difficult, as matrices in wavevector space be- 
come more complicated, and approximate diagonaliza- 
tion methods such as those employed in Refs. ( fOpn) 
cannot be used. The only existing microscopic numerical 
self-consistent calculationsrHri addressing the proximity 
effect at an F/S interface, are based on an extended Hub- 
bard model in real space. This is adequate only when the 
coherence riength is very short, and the material param- 
eters usedu were unrealistic.. Analytic work is similarly 
hampered. The traditionalcl way out is to conjecture a 
dependence of the proximity length on the exchange field, 
but the underlying assumption, while plausible, has never 
been proved and has been recently labelednJ as being just 
ad hoc. Physically, the spin imbalance in F results in a 
modified Andreev process, .sijQce the electron and hole 
occupy opposite spin bands .cJ The exchange field causes 
the quasiparticles comprising a singlet Cooper pair to 
have different wavevectors, so that the pair anmJitude 
in the ferromagnet becomes spatially modulated.E£l Such 
oscillations were first investigated long ago by Fuldc and 
Ferrellcj and Larkin and Ovchinnikov.tZl The resulting 
oscillations in F(r) induce oscillations (about the nor- 
mal state value) in the local density of states (DOS) as 
a function of distance from the interface. These oscil- 
lations have been studied theoreticallycSl (but non self- 
consistently) by using the Eilenberger equations, and 
good agreement was found with experiment .E3 The Us- 
adel equations revealed similar behavior Ej Quasiclassical 
approaches are clearly no substitute for a microscopic 
analysis which is required to study the case when the 
correlation length is not small. 

Thus, the theoretical situation for the F/S interface 
problem is unsatisfactory. In this paper, we attack this 
problem by obtaining numerical, fully self-consistent so- 



lutions for the continuum BdG equations for a ferromag- 
net in contact with an s-wave superconductor. Our nu- 
merical iterative methods overcome the technical difficul- 
ties associated with the different Fermi wavevectors, al- 
luded to above, and allow us to focus on the case of longer 
superconducting coherence lengths in the clean limit. We 
are able also to allow for different bancbsjidths in the two 
materials (Fermi wavevector mismatchEiJ) . The full BdG 
equations that are our starting point provide a rigorous, 
microscopic method for studying inhomogeneous super- 
conductors and their interfacial properties, and have the 
advantage that their solution provides the quasiparticle 
amplitudes and excitation energies. The resulting wave 
functions and energies are used to compute physically 
relevant quantities. We extract the relevant lengths from 
analysis of F(r) and investigate the local DOS as a func- 
tion of position on both sides of the F/S interface. Our 
results put the entire theory of the F/S proximity effect 
on firmer grounds, confirm some of the features previ- 
ously obtained approximately, and uncover new ones. 

This paper is organized as follows. In Sec. |l| we 
introduce the spin-dependent BdG equations, and the 
methods we employ to extract the pair potential, the pair 
amplitude, and the local DOS. In Sec. Ill we discuss the 



physical parameters we will use and present the results. 
Finally in Sec. 
future work. 
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II. METHOD 

In this section we present the basic equations we use 
for a system containing a ferromagnet/superconductor 
(F/S) interface and the methods we employ for their self- 
consistent solution. After self-consistency for the pair 
potential is achieved, we can then calculate other phys- 
ically relevant quantities such as the pair condensation 
amplitude and the local DOS. 

The system we consider is semi-infinite and uniform in 
the X, y directions and confined to the region < z < d, 
with the F/S interface located dX z = d' and the super- 
conductor in the region z > d' . We will take here d and d' 
larger than the other relevant lengths in the problem in 
order to study the interface between two bulk materials. 

We begin with a brief review of the starting equations 
in order to clarify our notation and conventions, including 
spin and choice of parameters. For a spatially inhomoge- 
neous system, a complete description of the quasiparticle 
excitation spectrum along with the auasiparticle ampli- 
tudes is given by the BdG equationaS. In the absence of 
an applied magnetic field, the system is described, using 
the usual second quantized form, by an effective mean 
field Hamiltonian, 

Heff = ^y rfV{^t(j.)Ho(r)^,(r) + V^i(r)/i,,,(r)^,-(r) 

a.a' 

+ i,7..,[A(r)^t(r)^t,(j.) + A*(r)V^.(r)^.,(r)]}, (1) 



where A(r) is the pair potential, to be calculated self- 
consistently, greek indices denote spin, rj = iay (the ct's 
are the usual Pauli matrices), and h{r) — —hoazQ{d'~z) 
is the magnetic exchange matrix. The step function in 
this term reflects the assumption that the exchange field 
arises from the electronic structure in the F side. The 
single-particle Hamiltonian is given by. 



Wo(r) = --^V2 + C/o(r)-A*. 
2m 



(2) 



Here fj, is the chemical potential, Uq is the spin indepen- 
dent mean field term, and we have set ?i = fc^ = 1. 

The BdG equations are derived by setting up the diag- 
onalization of the effective Hamiltonian via a Bogoliubov 
transformation, which in our notation is written as 

^T (i") = J2 [""T (r)7n - <T ('■)^n] ' (3a) 

n 

$1 (r) = J2 [""i (^)>' + <i ('■)^"] ' (3b) 



The Hamiltonian is translationally invariant in any 
plane parallel to the interface, therefore the component 
of the wavevector perpendicular to the z direction, k±, 
is a good quantum number. We can then write 



u„,(r) = <(z)e'''-■^ 
«„,(r)=<(z)e''^--. 



(7a) 
(7b) 



where k^ = (fcx, ky,0). Eqs.(0) then become one dimen- 
sional BdG equations. 



1 d^ 1 



-l-e^ + h^a{z) -/i <(2;) 



1 92 



(8a) 



(8b) 



where 7 and 7^ are Bogoliubov quasiparticle annihila- 
tion and creation operators respectively, and n labels 
the relevant quantum numbers. The quasiparticle am- 
plitudes Una and Vna are to be determined by eequir- 
ing that Eqs.(||) diagonahze Eq.(|l]). The resultingP BdG 
equations read, 

{Ho + K„)una + ^ P0-0-' A(r)w„cr' (r) = e„w„o-(i"), 

a' 

(4a) 
-(TYo + h„„)vna -|-^Pcrcr'A*(r)w„o-'(r) = e„w„o-(r), 

(4b) 

where p = ^x^ and the e„ are the quasiparticle en- 
ergy eigenvalues measured with respect to the chem- 
ical potential. Equations (0) must be supplemented 
by the self consistency condition for the pair potential, 
A(r) — g{r){'ip'^{r)-ipi{r)), which in terms of the quasipar- 
ticle amplitudes reads, 

^(r) = ^ E P..'E'"«-(^)<-' (^) tanh(e„/2r), (5) 

where g{r) is the effective superconducting coupling. We 
take this quantity to be a constant in the superconduc- 
tor, and to vanish outside of it. This is analogous to 
the assumption made for h. Our method does not re- 
quire that a small nonzero value of g be assumed in 
the non-superconducting side. The prime on the sum 
in ^ reflects that the sum is only over eigenstates with 
|en| < W£), where ujd is the cutoff (Debye) energy. The 
normalization condition for the quasiparticle amplitudes 
in our geometry is. 
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d\[\Una{r)\'' + \Vnair)\' 



1. 



(6) 



where ej_ is the transverse kinetic energy, and we have 
absorbed the mean field term by a shift in the zero of the 
energies. One can assume A(z) to be real without loss of 
generality. 

We can now solve Eq. (^) by expanding the quasipar- 
ticle amplitudes in terms of a complete set of functions 

(I^Tniz), 



K{z) 



N 

E 

771—1 

N 



Km(t>rrL{z), 



"Ziz) = E ^nm(t^"^iz)- 



(9a) 
(9b) 



A set of functions appropriate for our setup and geometry 
is that of the normalized free particle wavefunctions of a 
one-dimensional box, 



(j^miz) 



-sm{kmz), 



rmr 



If there was only one Fermi wavevector in the problem, 
the upper limit in the sum, N, would be-determined by 
that wavevector and ujd in the usual way.E3 But since this 
is not the case some care is required. The appropriate 
cutoff for this problem is given by 



N = [{kFxd/7r)y/l + LJD/fA, 



(11) 



where kpx is the largest Fermi wavevector in either the S 
or F side (see below) and the brackets denote the integer 
value of the expression they enclose. In a similar way, we 
can also expand the pair potential. 



N 






(12) 



After inserting these expansions into Eqs.(^ and making 
use of the orthogonahty of the chosen basis, we obtain the 
foUowing equations for the the matrix elements, 



[2^ + ^^1 ""9 " ^ [^^^'' ~ ^""^^««' + ^^' 



Oqq> 



a' p.p' 
r fr2 



L2m 



+ eAv: 



Iq + zZ [i^FM 



ha<7)Fqq' + EpsSq 



/ , / ^ Paa' ^pJpp'q '^np' ~ ^n^nq- 
a' p,p' 



"-nq' 

(13a) 

(7 

nq' 

(13b) 



In writing each term in Eq. (O) we have taken care 
to measure the chemical potential from the same origin 
(bottom of the same band) as the corresponding ener- 
gies. Because of the magnetic polarization and possible 
differences in carrier densities between the ferromagnet 
and superconductor, there are up to three different Fermi 
wavevectors involved in the problem, the two correspond- 
ing to spin up and and spin down on the F side, and 
one in the superconductor. On the F side we have in- 
troduced EpM through kp^/2m = Ep-^ = EpM + h^, 
kp,/2m = Epi = EpM — ho- On the S side, we have 
kpg/2m = Eps, where Eifs is the appropriate band- 
width. It has been shown,E2l that Fermi wavevector mis- 
match in F/S tunneling junction spectroscopy leads to 
nontrivial differences in the conductance spectrum. The 
matrix elements in (O) are given by, 



Fqq' — Cq'-q{d ) - Cq'+q{d ), 



(14a) 
(14b) 



Jpp'q 



E. 



q+p' 



-p{d) - Eq+pi^p{0) + Ep+.p>^q{d) 



Ep+p'-q{0) + Ep+q-p'{d) - Ep+g^p/{0) 



Ep+p'+q{d) + -Ep+p'+q(0) 



(14c) 



where we have defined Cm{z) = sin(fc,„z)/(7rTO), 
Em{z) = cos{kmz) / {nm) , for m ^ 0, and Eq{z) = 1. 
The self-consistency condition now reads, 

^« = f E ^PP'? E T!p--'<p<p' tanh(e„/2T), (15) 

p,p' a. a' n 

where the quantum numbers n include e_L and a longitu- 
dinal index m, the sum being limited by the restriction 
mentioned below Eq.(|^), and we have, 

Kppiq = -= Eq+pi^p{d) — Eq+pi^p{d ) + Ep+pi^q{d) 

\J id^ 

~ EpJ^pi -q\d ) + Ep^q^pi (d) — Ep^q_pi[d ) 

- Ep+pf+q{d) + Ep+pf+q{d') . (16) 



Finally, the normalization condition, Eq. 
the expansion coefficients, is 



in terms of 



EElKn 






^] = 1. 



(17) 



It is very difficult to solve Eqs.(|l3|) numerically as they 
stand, for large sizes. The required effort can be con- 
siderably reduced by solving for w-Ji„,u4„ only, allowing 
for both positive and negative energies. The solutions 



for u\ ,v\„ are then obtained via the transformation: 



U. 



T 

nq 



nq^ nq 



— > — e„. This simplifica- 
tion follows from the form of the exchange matrix, be- 
low Eq.(|^). Formally, the exchasge field breaks the ro- 
tational invariance in spin spaceJ23 however, there are no 
spin flip effects, so that the four equations (jl^) split into 
two equivalent sets of equations. 

For any fixed e± we can now cast Eqs. (^3|) as a 2N x 
2N matrix eigenvalue problem. 



H+ D 
D H- 



*n = e„*„. 



(18) 



where ^„ is the column vector corresponding to '^J-^ — 



«i 



' • ■ ■ I "nAT; 



w„]^, . . . , w„^). The matrix elements are 



H+q, = [2^ + £A-]5qq, - Ep^Fqq, - EpgSqq, , (19a) 

k"^ 

H-q, = -[^ + £i] V + EpiFqq, + EpsSqq, , (19b) 

pqq' ■ i^9CJ 



ii' ^ Z_^ ^pJj 



The basic method of self consistent solution of Eqs. 
( p8| ) and (Qq) works as follows: we first choose an initial 
trial form for the Ap. We then find, by numerical di- 
agonalization, all the eigenvectors and eigenvalues of the 
matrix in Eq. (Hq) , for every value of £±_ consistent with 
the energy cutoff (see Eq. (|ll|)). The formally continu- 
ous variable ej_ is discretized for numerical purposes. The 
calculated eigenvectors and eigenvalues are then summed 
according to Eq.(|l5|), and a new pair potential is found. 
This new pair potential is then substitutedo into the 
entire set of eigenvalue equations, and a new set of eigen- 
values and eigenvectors is obtained, from which in turn a 
new pair potential is constructed. The whole process is 
repeated until convergence is obtained, that is, until the 
maximum relative change in the pair potential between 
successive iterations is sufficiently small (see below). As 
an initial guess for the pair potential one can use, in the 
first instance, a step function of the bulk value, Aq, in 
the superconductor. The initial Ap are then obtained 



by inverting (12). After self-consistent results for Ap for 
one set of parameter values have been obtained, those re- 
sults can be used as the initial guess for a case involving 
a nearby set of parameter values. This process reduces 
the number of required iterations considerably. The final 
self-consistent result is insensitive to the initial choice. 
By using these methods, it is then possible, as we shall 



see, to obtain results even when the coherence length is 
long. 

This general procedure immediately yields the self con- 
sistent results for the pair potential. As mentioned in the 
introduction, this quantity gives valuable information re- 
garding superconducting correlations on the S side only, 
since it vanishes on the F side where g(r) = 0. Insight 
into the superconducting correlations on the F side, and 
the extraction of the proximity effect in thfi_ferromag- 
net, is most easily obtained by consideringpo the pair 
amplitude. 



F(r) = A(r)/g(r), 



(20) 



which has a finite value on both sides of the interface. 
One can also study proximity effects through another 
quantity which is directly related to observation. This 
is the local density of states (DOS), given byc3 

N{z, e) = ^ ^'[ML(z)'^(e - en) + vl,{z)Sie + €„)]• 



(21) 

In the next section we will first consider the relevant 
set of dimensionless parameters in the problem, and how 
we implement the general procedures discussed above for 
a wide range of the values of these parameters. We then 
discuss our results, and investigate the length scales rel- 
evant to the variation of the pair potential, the pair am- 
plitude and the DOS. 



III. RESULTS 

Before discussing our numerical techniques and results 
for the model outlined in the previous section, we have 
to introduce a convenient set of dimensionless parame- 
ters for the problem. First, there are two dimensionless 
ratios arising from the three material parameters Eps, 
EpM and /iq- We choose the ratio / = ho/EpM as the 
dimensionless exchange field parameter we will vary to 
study different degrees of polarization for the F side. / 
varies between / = when one has a normal (nonmag- 
netic) metal and 1=1, the half- metallic limit. In this 
work we choose the second ratio so that Ep^ /Eps = 1 at 
the value of / under consideration. Next, we have to con- 
sider the superconducting parameters. We have chosen 
to present here results for T = 0, postponing the study 
of temperature effects for future work. We then need to 
specify the dimensionless Debye frequency u = ud/ Eps 
and the dimensionless length scale kps£,Q, where ^o is the 
usual zero-temperature coherence length related to other 
quantities by the BCS relation, kps^o = (2/7r)(i?Fs/Ao). 
Throughout, we will keep the relatively unimportant pa- 
rameter Lo fixed at 0.1, and present results for two differ- 
ent values of kps^o, 50 and 200. Thus, our method can 
handle coherence lengths two orders of magnitude larger 



than what has been achieved through the usecJ of tight 
binding methods. 

We also have to consider the purely computational pa- 
rameters. These are determined by the overall size of the 
system, measured in terms of kpsd, and the ratio d! / d. 
Our two choices of ^o demand different system sizes, since 
the length scale over which the pair potential reaches its 
bulk value in S is determined (see below) by ^q- Thus, we 
need d 3> Co in order to study an interface between bulk 
systems. Thus, we take kpsd = 1000, d'/d = 600/1000, 
for kpsS,o = 50, and kpsd = 1700, d'/d = 750/1700 for 
kps£,o = 200. These slab widths allow us to investigate 
fully the bulk proximity effects that occur on both sides 
of the interface. 

The computational work required is chiefly determined 
by the system size. As outlined in the previous section, 
we must numerically diagonalize the Hamiltonian matrix 
and calculate the eigenenergies and eigenvectors for each 
e^. Each value of e± requires diagonalizing a matrix 
of size 2N x 2N, where N is defined in Eq. (|ll|). For 
the large values of d required by our assumed values of 
Co this matrix size exceeds 1100. The number of dis- 
cretized transverse energies, N±, must be chosen large 
enough so that the results are not affected by it. The re- 
quired value depends on the quantity being studied. For 
A(z), and F{z), a value of N± = 500 was found to be 
sufhcient even for the longer coherence length. For the 
local DOS, we used iVj^ = 1000 in both cases. These 
diagonalizations must be performed at each step in the 
iteration process described below Eq. (|l9[) . The basic di- 
agonalization process employed a procedure whereby the 
symmetric matrix, Eq.(|lq), is transformed into tridiago- 
nal form, and then the eigenvalues and eigenvectors are 
computed by the QLE3 algorithm. The iteration process 
was concluded when the maximum relative error between 
successive iterations of the pair potential at any point 
was less than 10"'^. A smaller relative error would re- 
quire more computation time, but we verified that no 
appreciable difference in the results ensued. A number 
of checks were performed, including reproducing the cor- 
rect wavefunctions and energies for the limiting case of a 
single semi-infinite superconductor, ferromagnet or nor- 
mal metal, and also verifying that in the limit of an 
entirely supepionducting sample the correct finite size 
oscillationsEHH of the pairj^otential were obtained, with 
the correct Co dependences 



A. Pair potential 

We begin by presenting in Fig.n^ our self consistent re- 
sults for the pair potential, A(z) (normalized to the bulk 
value Ao), which we plot as a function of the dimension- 
less variable Z = kpsz. In the four panels of the left 
column we show results for kps^o = 50 for four evenly 
spaced values of / ranging from zero to unity. In the 
corresponding panels in the right column we have results 
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FIG. 1. The self consistent pair potential A(z), normalized 
to the bulk value Ao, is plotted as a function of dimensionless 
distance Z = kpsz. The left column is for fcF^o ~ 50, with the 
F/S interface at Z = 600. The right column is for kp^o ~ 200, 
the interface is then ai Z = 750. In both cases results for the 
same four exchange fields (indicated by the labels) are shown. 



FIG. 2. The pair amplitude F{z) (defined in Eq.(g|)), nor- 
malized to its bulk value in the superconductor, plotted as a 
function of dimensionless distance Z = kpsZ- Results are for 
the same coherence lengths and exchange fields as in Fig.g], 
and with the same panel arrangements. 



for kFsS.0 = 200 at the same values of the exchange field. 
The pair potential always vanishes on the F side since 
we assumed ^(r) = in that region. All of the panels 
show that on the S side, the normalized pair potential 
rises near the interface and then eventually reaches its 
bulk value over a length scale determined by the coher- 
ence length ^0- Comparing the top panels in each col- 
umn, where / = 0, with the others in the figure, where 
the exchange field can be large, we see that for all four 
values of / and a given value of ^Oj the characteristic de- 
pletion near the interface is nearly independent of /. It 
can be concluded therefore, that the magnitude of the 
exchange field has little effect on A(z) and that the ef- 
fective coherence length in the superconducting side of 
the F/S interface is only an extremely weak function of 
the strength of the ferromagnetic exchange field. Similar 
findings were obtained in Ref. 21, in the short ^o limit. 



The general shape of the curves is the same for both 
values of ^O: indicating that the effect of this quantity 
is merely a rescaling of the relevant length which gov- 
erns the interface depletion. Near the surface-vacuum 
boundary in S, the pair potential exhibits atomic scale 
(a distance of order l/kps) oscillations as seen in previ- 



ous worko, as a result of pair-breaking by the surface. 



B. Pair amplitude 

The above study of A (2:) illustrates the detail, and 
quality of the results. However, since A(z) vanishes in 
the F side, this quantity cannot be used to study su- 
perconducting proximity effects in the magnet. For this 
purpose we now turn our attention to-the pair amplitude 
F(z), a quantity that directly reflect^ the superconduct- 
ing correlations in both F and S. The main panels in Fig- 
ure]^, which repeat the arrangement of Fig. y| show eight 
sets of results for F{z), four for each of our two values 
of kps^o, for the same values of / as in Fig.|l]. We have 
normalized F{z) to its bulk value in the superconductor. 
In the S region the curves are the same as those for the 
corresponding A(z), seen above in Fig. 111. 

Turningour attention to the F side of the interface (z < 
d') in Fig.||, we first look at the normal metal Ikpt (/ = 
0) in the top panels. We see that as expectedclEa F{z) 
decays only extremely slowly in the normal region at T = 
0. In effect, there is no mechanism to dismpt the Cooper 
pairs from drifting across the interfaceEj, therefore the 
decay is very slow and occurs over a length scale that 
is much larger than ^o- The most rapid change occurs 
near the interface, where F{z) decays very quickly before 
flattening out. 

In the remaining panels of FigJ3 the effects of a finite 
exchange field are seen. The situation is now very dif- 
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FIG. 3. Detail of the behavior of the normahzed pair am- 
phtude F{z) near the interface, on the magnetic side. The 
six panels shown correspond to the lower six panels of Fig. 
0, but the horizontal scale is expanded so that the oscillatory 
behavior can be seen. 



ferent and F{z) decays to zero rather quickly close to 
the interface, with a slope that increases with larger /. 
We will see below that the length that characterizes this 
fast decay varies approximately as 1//. This 1// be- 
havior was suggested long agoQ on the intuitive grounds 
that the exchange potentials seen by up and down spin 
quasiparticles differ by ±7, but this argument has been 
criticizedcj as being merely an ad hoc assumption. Our 
results show that the intuitive assumption gives the cor- 
rect result. However, this fast decay is far from the whole 
story, as slightly away from the interface a much slower 
oscillatory behavior can be seen (note in particular the 
/ = 1/3 panel). This is not a finite size effect. We have 
replotted this behavior in an expanded horizontal scale 
in Fig.^ The wavelength of these oscillations clearly de- 
creases with increased /. Furthermore, the magnitude 
of F{z) also attenuates with increasing /. This is in 
qualitatiije agreement wittt-past works employing tight- 
bindingtil or quasiclassicafia methods. 

Before we consider this behavior in more detail let us 
examine the ^o dependence of these results for F{z), by 
comparing the right and left column of Fig. 0. The 
spatial extent in which the changes in F{z) take place 
is greater in the right column, since we are dealing now 
with a length scale given by a longer ^o- Apart from that, 
the differences are hardly discernible, the only exception 
being the very slow decay for / = 0, where the difference 
can be attributed to the smaller value of d'/'Co ~ 4 com- 
pared with d'/^o ~ 12 for the case of kps^o = 50. Thus 
we conclude that the role of ^o is, in this range, that of 
setting an overall scale. This should hold only when i^o is 
much larger than the microscopic lengths in the problem 
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FIG. 4. Example of a fit of the results for the normalized 
F{z) (solid curve) to the expression in Eq. (p3) (dashed line). 
The data displayed are for kps^o ~ 50 and / — 1/3 (as shown 
in Figs. and H 



and smaller than the geometrical dimensions. It should 
break down in any other case. The exchange field tends to 
disrupt superconducting correlations over a length scale 
that is typically much smaller than ^q, so that the oscil- 
lations and characteristic decay of F(z) in the magnetic 
region are nearly independent of the ^o considered here. 
We are then led to conclude that when there is an ex- 
change field present, there are two phenomena to consider 
in describing the spatial variations of F{z) in the ferro- 
magnet. The first is the short distance decay at the inter- 
face, to the point at which the pair amplitude first goes to 
zero. This is the region where F{z) changes most rapidly. 
This decay can be characterized by a length scale which 
we will denote by ^i , and define as ^i = d' — zi, where zi 
is the first point inside the magnet where F{zi) = 0. The 
other important phenomenon is the damped oscillations 
of F{z) in the region z < zi (Fig. 0). These oscillations 
cannot be fit to an exponentially damped form. Instead, 
we find that in all cases a much better fit to our results 
is afforded by the following expression: 



F{z) = a 



sin[(z ~ dO/6 



(22) 



where a is a constant, and the characteristic length ^2, 
which in principle must be distinguished from ^i , can be 
extracted from the results. Since the previously defined 
length, ^1, is small, the expression ( p2| ) is valid for most 
of the ferromagnet region. To illustrate the range of its 
validity, in Fig.^ we give one example of a fit of the form 
Eq.(p2[) to the pair amplitude. We see that Eq.(p^) is an 
adequate fit for the oscillatory region, however, within 
a distance ^i of the interface, Eq.(p2|) breaks down. At 
this point, F{z) rises upwards monotonically to match 
its value at the interface. In the spatial region where 
Eq.(p2|) is valid, the quality of the fits deteriorates some- 



what for larger exchange fields (0.4 < / < 1) because 
the spatial modulation of F{z) slightly deviates from the 
simple periodic sine curve given by Eq.(p^. This small 
discrepancy can be glimpsed in the lower panels of FigJ^. 
The spatial structure becomes slightly nonperiodic, but 
overall the functional form given by Eq. (|2^) is still satis- 
factory. 

The oscillatory behavior of the-pair amplitude as given 
by (pa) is physically the resultcJ of the exchange field, 
which creates electron and hole excitations in oppo- 
site spin bands. The pair amplitude involves products 
of these particle and hole quasiparticle amplitudes (see 
Eq.(p|)). The superposition of these wavefunctions then 
creates oscillations on a length scale set by the difference 
between the spin up and spin dpwn wavevectors in the 
ferromagnet. One then expectsE3 a decay of the form 
( P2| ) with ^2 ~ (^FT ~ kpi)^^- We can then write: 
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where in the last step we have used, as previously men- 
tioned, Epi/Eps = 1- For small /, Eq.(2^) can be sim- 
plified to kFs£,2 ~ 1/^, showing that 1^2 is then inversely 
proportional to the exchange field. At larger values of 
/ there are deviations, but these are small since in the 
7=1 limit both Eq. (|2^) and the approximate expres- 
sion coincide. These oscillations are also related to those 
responsible for oscillatory coupling in structures involv- 
ing magnetic layers and superconducting spacers J23 and 
the nonmonotonic behavior in the critical tempersiture, 
Tc, versus the F-layer thickness in S/F/S junctions.L3 In 
particular, the sign change in the pair amplitude has the 
same physical origi n a s t he so called "tt phase" that exists 
in F/S multilayers ,E3cil and the nonmonotoniiL variation 
of the Josephson current with exchange ficld.Cj 

Having introduced the two length scales ^1 and ^2 char- 
acterizing the superconducting proximity effect in the 
magnetic region, it is useful to compare their magnitude 
and behavior as functions of I. The result of doing this 
is shown in Fig.0. Data at additional values of /, not 
displayed in previous Figures, is included. For compari- 
son, Eq.(g^) is shown as the solid curve. We find that ^2 
follows very closely the expected theoretical expression, 
and that the other length, ^i(/) ~ S,2{I\i This is be- 
cause, as mentioned above, the expressionu kps^i = l/I 
nearly coincides numerically with the more complicated 
result for ^2 as given above. Thus it turns out that the 
fast decay and the spatial period of the oscillations are 
characterized by lengths that are virtually identical. 



C. Local density of states 

To further investigate the F/S proximity effects, we 
focus now on another experimentally accessible quantity. 



LL 




FIG. 5. Exchange field, /, dependence of the lengths 
^i,i = 1,2, defined in the text. The circles are ^i, and the 
squares represent ^2. The curve is the expression in Eq.(E3|). 
The results plotted are for kps^o ~ 50 but these quantities 
are nearly independent of go • 



the local DOS. Advances in STM technologya have made 
it possible to perform localized spectroscopic measure- 
ments with atomic scale resolution. We therefore present 
now the local DOS as a function of energy and position, 
as calculated from Eq.([2l|) and the self-consistent spec- 
tra. All results below are normalized to the normal-state 
DOS in the S side and convolved with a Gaussian of width 
O.OIAq, to eliminate the spectrum discretization result- 
ing from the finite size of the computational sample. We 
focus only on results for positive energies, since those for 
negative ones can be obtained by symmetry. We plot the 
results in terms of the normalized energy variable e/Ag. 
The locations chosen are given by the dimensionless po- 
sition measured relative to the interface in terms of the 
quantity Z' = kps{z — d'). Thus, a positive value of Z' 
denotes a location within the superconductor. 

We consider first the limit where the exchange field / 
is zero. In Fig.H, we show the DOS for four different posi- 
tions at each of the two values kps^a = 50 (left column) 
and kps£,o — 200 (right column). The three top rows 
in Fig.H show the DOS on the S side. For the shorter 
coherence length results for the locations Z' — 50, 100, 
and 200 are shown. These are multiples of the coherence 
length, and the same multiples are shown in the right 
column. Several pronounced peaks are visible inside the 
gap, due to a finite number of bound states existing for 
e/Ao < 1. These states were predicted long ago in a 
non self-consistent treatment by de Gennes and Saint- 
Jamesllil. These peaks diminish at greater distances in- 
side the superconductor. On the corresponding panels in 
the right column, we see that the number of de Gennes 
Saint- James peaks have been reduced. This is because 
the number of bound states depends upon the coherence 
length ^Oj aSj-well as on the superconductor and normal 
metal widthst3. In general, the number of such peaks de- 
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FIG. 6. Normalized local DOS (see Eq.(|2l|)), plotted versus 
the dimensionless energy e/Ao at / = for kps^o = 50 (left 
column) and ^fsCo ~ 200 (right column). The values of the 
dimensionless quantity Z' = kpsiz ~ d') shown by the labels 
are positions relative to the interface. Each position displayed 
is a multiple of the coherence length: from top to bottom 
the rows correspond to Z' = ^o, Z' = 2^o, Z' = 4^o, and 
Z' = -2Co. 



creases as £,o/d' increases. The patterns seen at e/Ap > 1 
are discussed below. 

On the normal metal side, we see on the bottom panels 
of Fig.^, that there is no evidence of a gap, but a pattern 
of jagged peaks appears in the DOS for e/Ag ^ 1. At 
larger energies, interference patterns are seen, similar to 
those in the S side. At longer coherence lengths this 
pattern is more coarse. This coarseness (which is also 
seen in subsequent figures) arises from the finite value of 
N±. If this quantity is increased, the pattern becomes 
smoother and more regular, as in the left column. The 
remaining regular oscillations ultimately vanish as d and 
d' tend to infinity. We have chosen to display only one 
position in the F side for J = since the overall behavior 
is nearly identical for all points in the normal metal. This 
is in agreement with our observation in connection with 
the / = panels in Fig.0, that the pair amplitude has a 
very slow rate of change. 

We now turn to the case of a finite exchange field. Fig- 
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FIG. 7. Normalized local DOS for / = 1/3 at four positions 
in the ferromagnetic side, near the interface. The left column 
corresponds to kFs£,o ~ 50 and the right one to kps^o ~ 200. 
Z' is defined as in the previous Figure. 

ure 1^ shows the DOS for kFs£,o — 50 (left column), and 
kFs£,o = 200 (right column) at / = 1/3 for four posi- 
tions very near the interface, within the magnetic mate- 
rial. This is done to illustrate how changes in the local 
DOS with distance are correlated with the rapid change 
in F(z) near the interface. Consider first the distance 
Z' = — 5 (top panels). This corresponds to the location 
where F{z) has its more prominent minimum (see Fig. 
g). There is a weak minimum for the DOS at e/Ag = 0, 
which is more prominent at the smaller ^Oi and with in- 
creasing energy the DOS rises, until about e/Ag « 1 at 
which point a peak occurs. For energies larger than Aq 
the DOS quickly settles down to its normal state value, 
unity in our normalization. At Z' = —4, as F{z) be- 
gins to rise, we see, focusing on the range of energies 
less than Aq, that the minimum of the DOS has begun 
to shift away from zero. At Z' = —3, in the next row of 
panels, the DOS has now a marked minimum at finite en- 
ergies within the gap, at e/Ag ~ 0.6. The next position 
(last row) in Fig.|^ shows a clear minimum of the DOS at 
energies just below the gap. By comparing the top and 
bottom rows of Fig.R, we see that (for e/Ag < 1) what 
were once dips and peaks in the DOS have now reversed 
roles. Fig.H shows that the length characterizing the fast 
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FIG. 8. Local DOS at three positions inside the supercon- 
ductor, for / = 1/3 and kFs£,o = 50. The curves shown 
correspond to (from top to bottom at small energy) Z' = 50, 
Z' = 100, and Z' = 200. These are the same values as in the 
top three panels of the left column of Fig. o. 



rise of F{z) is kps^i ~ 3.5 at / = 1/3. The DOS starts 
the reversal process, as the interface is approached, at 
around Z' « —3.5, as seen in Fig. M. The similarity be- 
tween right and left columns in this Figure reflects that 
the length scale ^i, defining the inversion point, is the 
same in both cases. The behavior of the DOS at larger 
values of \Z'\ is qualitatively similar, but as the oscilla- 
tions die down it becomes much less discernible. 

We show also (see FigJ|) the local DOS for three dif- 
ferent positions in the superconductor side for the same 
7 = 1/3 and the shorter ^o- The de Gennes-Saint James 
peaks are now gone, with just a hint of small structure 
for s/Aq < 1 remaining at Z' — 50. This structure starts 
to become washed out at a distance of about 2^o from the 
interface. Finally, at Z' — 200, the DOS is of the familiar 
BCS form, with a well defined gap and pronounced peak 
at s/Aq — 1. In what follows, we focus only on the fer- 
romagnetic region, since the overall behavior of the DOS 
in S at larger values of / is quite similar to that seen in 
Fig.|. 

In Fig.|, we show the DOS for / = 2/3. As in Figj^, 
we consider four spatial positions for each value of ^o, 
however the range is now closer to the interface, since 
the larger exchange field reduces the spatial extent of 
the superconducting correlations and the length ^i. Be- 
ginning at Z' — —4, Fig.^ (top) illustrates the formation 
of a small dip at low energies, and a continual rise up to 
s/Aq « 1, after which we recover the bulk DOS limit for 
a ferromagnet with this polarization. With our normal- 
ization, this value is smaller than unity. This is due to the 
decrease in the number of spin down states with increas- 
ing exchange field. In Fig.0 (second row), the minimum 
in the DOS has moved, while the peak still remains at 
s/Ao « 1. At Z' = -2, (third row) the DOS is aheady 
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FIG. 9. Normahzed local DOS for / = 2/3. The panel 
arrangement is the same as in Fig. m. 



rising upwards at low energies. The previous dip in the 
DOS has shifted to a higher energy, while a peak forms 
around zero energy. We again find consistency with the 
fi values given in Fig.^, where for / = 2/3, kps£.i ~ 2.2. 
Thus, as in the previous case, a reversal of the DOS be- 
havior occurs in the ^i range. The bottom panel of Fig. ^ 
shows the DOS at Z' = —1. We see that the zero energy 
maximum has increased slightly from the previous row 
and the minimum has shifted to energy Aq. Again, this 
qualitative behavior is independent of ^o reflecting the 
independence of ^i from ^o- Thus, we see here the same 
behavior we found for / ~ 1/3 the only change being the 
different value of ^i . 

We flnally consider in Fig.nn the DOS for a fully polar- 
ized (half metallic) ferromagnet (/ = 1.0). The locations 
Z' are the same as in Fig.pl The structure of the DOS 
at energies below the gap for all positions has become 
smoother. Because of the large exchange flcld, the rever- 
sal of the occupation of states occurs over a length scale 
^1 which is now small (see Fig.0). Based on the previous 
flt in Fig||, we find this point to be Z' « 1.7. Again, 
we find consistency between the pair amplitude and the 
DOS. Note that as one moves away from the interface, 
the DOS tends to 1/2 at higher energies. This is due to 
the total absence of the down band in this half metallic 
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FIG. 10. Normalized local DOS for a fully polarized ferro- 
magnet (7 = 1.0). Again, results for AfsCo ~ 50 are in the 
left column, and those for kps^o = 200 in the right column, 
and positions with respect to the interface are indicated. 



limit. These remarks apply to both values of ^o- 

The above pertains to the proximity effect in the F 
side. We now investigate further whether the presence 
of the ferromagnet has any effect on the superconduct- 
ing correlations. That any effect is small is shown al- 
ready by Fig. 0. The influence of the magnet on S 
will be reflected in a nonzero value of the difference in 
the density of states for spin-up and spin-down electrons, 
6N = Ni~Ni, where N-^ and N^ are the spin up and spin 
down terms in Eq. (E^) respectively. One might view this 
as a self consistent determination of an effective parame- 
ter I(z) which may extend into the superconductor. We 
focus here only on the case of a half metallic ferromagnet 
(/ = 1), and for illustration take kps^o = 50. Figure pTJ 
shows that there is in fact a small proximity effect into 
the superconductor, since very close to the interface the 
effective polarization is nonzero. This effect is however 
short ranged, and we see that it nearly dies out before 
Z' = 5. At very small exchange fields (of order of the 
superconducting gap) , we have found also a longer range 
proximity effect in the jSuperconductor, similar to that 
found for dirty systems.CJ 
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FIG. 11. Leakage of magnetism into the superconductor; 
the quantity plotted difference between spin up and spin down 
values of the local DOS, 5N = N^ — Ni, normalized as in pre- 
vious Figures. This quantity is displayed at several positions 
just inside the superconductor. All results shown in this Fig- 
ure are for kps^o ~ 50. The values of Z' are indicated in the 
panels. 



IV. CONCLUSIONS 

We have introduced in this paper numerical techniques 
to accurately and self-consistently solve the continuum 
BdG equations. We have shown how one can use these 
methods to perform a detailed study of F/S interfacial 
properties. Our procedures allow us to consider super- 
conducting and magnetic proximity effects in a bulk sys- 
tem containing an F/S interface, even when the super- 
conducting coherence length is orders of magnitude larger 
than the interparticle distance. In this work, we have 
used these techniques to investigate the proximity effects 
for a clean F/S system. We have extracted the relevant 
characteristic lengths through a careful analysis of the 
pair potential and the pair amplitude, and we have shown 
how, near the interface, the behavior of the pair ampli- 
tude correlates with that of the local DOS. Our work 
extends well beyond previous numerical computations in 
the tight-binding case, valid only for very short values of 
^0, and beyond theoretical work limited by quasiclassical 
approximations or restricted to regimes where the mean 
free path or ^o are very short. 

On the S side we have found, near the interface, a 
depletion of both the pair potential and the pair ampli- 
tude. This depletion extends over a length scale deter- 
mined essentially by ^o i and hence nearly independent of 
the exchange field /. For the bulk heterostructures we 
considered, the effect of varying ^o in the range studied 
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{kFs£,o = 50 to kpsio = 200) was an effective rescaling of 
the characteristic length that determines this depletion. 
In the F region, for finite values of the exchange field, the 
pair amplitude exhibited a sharp monotonic decline near 
the interface, followed by damped oscillations. The fast 
decay was found to take place over a length scale <^i ap- 
proximately inversely proportional to /, independent of 
the ^0, according the the expressiorjj kpsS,! ~ !//■ The 
oscillatory part of the spatial variation of the pair ampli- 
tude could be fit to a simple sine function with an ampli- 
tude decaying as the inverse of distance from the inter- 
face. We found that the spatial period of the oscillations 
is determined by the length difference ^2 = [kpi — kpi)^^ 
(the inverse of the difference between spin up and spin 
down Fermi wavevectors) provided that / is not too large. 
This is in reasonable agreement with previous theoretical 
expectationscj. We have presented extensive results for 
the local DOS, as a function of position and energy, as 
obtained via the self-consistent quasiparticle amplitudes 
and energies. The periodic sign change in the pair am- 
plitude is found to be correlated with oscillations in the 
local DOS relative to its normal state values. Finally, 
we verified also from the local DOS that the effect of 
the exchange field on superconducting correlations in S 
is minimal (although nonzero): the difference in the lo- 
cal DOS of spin up and spin down quasiparticles vanishes 
except very close to the interface, at least for / > 1/3. 

Clearly, the powerful methods and techniques for the 
self-consistent solution of the BdG equations presented 
here open new vistas and possibilities for use in the study 
of many other aspects of the F/S interface and similar 
problems. A thorough investigation of the physical quan- 
tities and characteristic lengths studied in this paper, in- 
corporating other parameter regimes and the effects of 
finite temperature is needed, and it can be straightfor- 
wardly carried out. Interface scattering, and supercon- 
ductors with nodes in the pair potential (unconventional 
pair potentials), can be also easily considered. Spin-flip 
effects, and disorder in both F and S materials can also be 
incorporated. By suitably changing the boundary con- 
ditions, self consistent solutions of the tunneling spec- 
troscopy problem in the long ^0 regime will be obtain- 
able. Our numerical methods are particularly suitable to 
the study of mesoscopic structures involving F/S multi- 
layers of differing thickness, where size effects may come 
into play. The study of tunneling phenomena in non- 
equilibrium situations is also feasible by extension of our 
method to the time-dependent BdG equations. 
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